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1. Introduction 

The wave propagation in bent waveguides has a long and rich history of research that 
dates back to Lord Rayleigh [1] and continues to the present days. Initially bends have 
been investigated in the framework of the electromagnetic theory, but more recently 
also quantum mechanical aspects attracted a lot of attention. The bends are popular 
subject of investigation, because they are typical elements incorporated into designs of 
waveguides. The computation of their properties to sufficiently high precision seems to 
be a difficult problem in the regimes of high energies and high curvatures even today. 

We are discussing a bend as a scatterer of non-relativistic quantum waves on a two- 
dimensional ideal straight waveguide as shown in figure 1. Such a structure is referred 
to as an open billiard. The past research of quantum aspects of bends can be separated 
into two branches. These are studies of bound states, their existence [2, 3, 4] and spectra 
[-] and the scattering properties, which are both reviewed in reference [6]. In order to 
describe quantum phenomena over our open billiard several approaches have been used 
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in the past: Green function approach [ ], finite difference mesh calculations and mode- 
matching techniques (MMT) using natural modes i.e. eigenfunctions of the Laplacian in 
the bend [9, 10, 11, 12, 5, 13] and other bases [14]. The work [ ] is particularly interesting 
as it raises the question on how to stabilise the calculations and gives a MMT method 
that is stable, but unfortunately a bit ambiguous. The MMT based on natural modes is 
called the modal approach and is the main topic of discussion in the present paper. The 
modal approach looks the most promising to deal with bends, because of its simplicity, 
power of interpretation and precision of results. But it also hides some problems that 
we examine here in detail. 




C 



Figure 1. A schematic picture of a finite bend of the inner radius r = q and the outer 
radius r = 1 on a straight waveguide of width a = 1 — q. 

Let us introduce the modal approach in our open billiard composed of a bend, 
with the inner radius r = q and the outer radius r = 1, and a straight waveguide of 
width a = 1 — g as shown in figure 1. The area of the billiard, denoted by can be 
separated into three sections: A - the left lead, B the bend and C - the right lead. We 
are searching for the wave function '0(r G f^), which solves the stationary Schrodinger 
(Helmholtz) equation on Q with Dirichlet boundary conditions reading 

- AV(r) = e^jiv) , i^lredu = , (1) 

where E = is the energy and k the corresponding wavenumber. The Helmholtz 
equation (1) is written in Cartesian coordinates r = {x^y) in the asymptotic regions A 
and C as 
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and in polar coordinates r = (r, (f>) across the bend B as 
52 Id 1 52 



9^2 ' rdr ' r2«-^2 /^ = ^ V-, ^\r=,,i = . (3) 

The modal approach suggests that we first solve Helmholtz equation (1) on each region of 
the open billiard separately. Thereby we obtain partial solutions called mode functions^ 
which are then used to describe the solution across the whole billiard. In the infinite 
straight waveguide and in the bend the mode functions are given by the ansatze 
ip oc u{y) exp{igx) and ip oc C/(r) exp(ii/0), respectively. The ansatze in equations (2) 
and (3) give the following equations for the mode functions: 

— + - g^)u = , u\y=o,a = , ?/ G [0, a] , (4) 

^ + ^^+(fc^-^)f/ = 0, t/U, = 0, rG[,,l], (5) 

where a = 1 — g is the channel's width and the scalars /a G C are called mode numbers. 
The equations (4) and (5) have a discrete set of solutions i.e. the mode functions and 
the corresponding mode numbers, which we refer to as modes. 

The modes in the straight waveguide and in the bend are denoted by pairs 
{9niUn{y)) and {iyp^Up{r))^ respectively, where n G N. The modes in the straight 
waveguide are explicitly written as 

= gl+ (^^) , Un{x) = y^sin (^^^) , (6) 

whereas modes in the bend are more complicated. They are discussed in Section 2, where 
we also show that mode numbers are either real or imaginary. The mode functions Un{y) 
and Up{r) are then used in the bases of functions in which we expand waves over parts 
of the waveguide. The basis in the straight waveguide is given by 
±/ X / xexp(±i^^x) 

en (r) = Un{y) — , (7) 

V yn 



and in the bend by 

/*(r) = UArf-!^^^ , (8) 



/Up 



where the sign ± labels the two directions of phase (probability) flux propagation. 
We define the square-root of a complex number z = |z| exp(i(/)), (j) G [0, 27r) as 
^/z = exp(i(/)/2). The basis functions are called wave modes or modes of the 

Laplacian. We distinguish two types of wave modes. The wave modes corresponding 
to real and imaginary mode numbers are called open modes or travelling waves and 
closed modes or decaying (evanescent) waves^ respectively. The wave-function '0(r) (1) 
in the entire open billiard region Q is expressed in terms of the wave modes as 

^(r) = 5^a+e+(r) + a;^e;^(r) , r G f^A , (9) 

n 

^(r) = yx+f+{r) + A-/-(r) , r G , (10) 
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t/.(r) = ^6+e+(r) + 6;e;(r), r G fic , (11) 

n 

where Qa,b,c are regions corresponding to sections A, B and C, respectively. The 
expansion coefficients a^, and are determined by the condition that the wave 
function ip{r) is smooth everywhere in $7, in particular on the boundaries between 
different regions Qa,b,c- The solution of the presented problem will be discussed in 
Section 3. 

The paper is organised as follows. In Section 2 we present a detailed study of the 
mode structure in the bend, which is closely related to the work of Cochran [15, 16, 17]. 
In comparison to the work of others, ours is directed more towards the application of 
the mode structure to scattering calculations. In addition, we write explicit formulae for 
the mode functions in the bend, where we give special attention to the closed modes. In 
Section 3 we outline a numerically stable MMT for calculation of the scattering matrix 
[18] of a single bend. The section 3 is concluded with the presentation of numerical results 
obtained by our method and compared to analytic estimates of the quantum transport 
properties of the bend. By considering the analogy between the quantum theory and EM 
theory we can connect our work to the EM wave propagation of longitudinal magnetic 
waves [9]. 



2. The cross-product of Bessel functions 

In this section we analyse the properties of the mode numbers and the corresponding 
mode functions for a given wavenumber k and inner radius q. The mode functions in 
the bend Up{r) are proportional to well known cross-products of Bessel functions [15] of 
the first kind, t J^^, and Bessel functions of the second kind, y^/, [ ] written as 

^.(fc, r) = J,{kr)Y,{k) - Y,{kr)J,{k) , (12) 



or 



sin(z/7r) 

where the allowed values of mode numbers v are determined by the Dirichlet 
boundary conditions 2'^,/c(?) = 0. In equation (13) we have used the relation Yy{z) = 
{Jiy{z) cos(z/7r) — J_^(z))/ sin(z/7r) valid for orders u ^ Z. The understanding of the mode 
structure in the bend is essential for calculations of the scattering over our open billiard 
in the modal approach. 



2.1. The properties of mode numbers 

The set of mode numbers at a given wave-number G M, > and inner radius 
q G (0, 1) is denoted by Aik,q = {^^ G C : Zy^^io) = 0}. The functions Z^^ki"^) are even 
Z_jy^k{r) = Ziy^ki"^) analytic in the order u [16]. These properties yield the following 
symmetry of the set of mode numbers: 

Mk,q = -Mk,q , = Mk,q • (14) 
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In addition we conclude that mode numbers are either purely real or purely imaginary 

Mk^q C MUiM. (15) 

The number of real modes is finite, whereas the number of imaginary modes is infinite 
at a finite wavenumber k. The proof of the later is given in Appendix A. The properties 
(14) and (15) enable a decomposition of M.k,q into two disjoint subsets of mode numbers 
laying on the positive M.k,q,+ and the negative M.k,q- real and imaginary axes: 

A^M,+ = ^ Mk^q : > or $5^/ > 0} , (16) 
Mk^q^- = {z/ G Mk^q : < or < 0} C -A^fc,^,+ . (17) 

It easy to see that Mk^q = U Mk,q-- We caU ?R:{Mk,q,^} and Q{Mk,q,^} the 

set of real modes and imaginary modes^ respectively. The number of real modes in the 
bend A^b = card 3?{A^/e,g,+} is equal or one more than the number of real modes in the 
straight waveguide Ag = [/ca/vrj: 



< Ab - A, < 1 



(18) 



where [xj denotes the largest integer smaller than x. Taking into account the analyticity 
of Zjy^ki^) in the order u and in the wavenumber k [16] we find that Ab can be computed 
in the semi-classical limit, for g 7^ 0, as 



A^b 



a 



ka 

7T Snqk 



+ 0{k-') 



(19) 



and Ag is asymptotically, as oc, close to Ab. The expressions (18) and (19) are 
explained in Appendix B. The asymptotic form of Zj^^^i"^) in the order parameter [20] 
reads as 



1 



Zr^Ar) = — [r^O(l) - r-^O(l)] . » 1 . 



(20) 



From equation (20) we learn that Zy^q) diverges exponentially with increasing order 
parameter on the real axis as 0(g~l^l) and oscillates along the imaginary axis. This 
bounds real mode numbers 3?{A^/e,r,+} from above and indicates that there is an infinite 
number of almost periodic imaginary mode numbers. 

The mode numbers v G M.k,r and consequently the mode functions can be in 
general obtained only numerically. We use different approximations of mode numbers 
to improve their numerical computation. By using the Debye approximation of Bessel 
functions for imaginary orders = iy (^ G M) valid for + {qk^ ^ 1 and the Dirichlet 
condition Z^^q) = we obtain the following relation 



Vy'^ + - Vy^^-Jkqy + y log 



y+ ^/y^ + k^ 



Trn , n G N 



(21) 



y + ^F+TM^. 

The solution of the equation (21) in variable y ^ k represents an asymptotic 
approximation of imaginary mode numbers Pn = Wn and is written as 

{kaf 



\ogq\ 



Aim 



+ 0(n-^) , n > 1 . 



(22) 
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The first term in equation (22) is already well known and can be also obtained from 
equation (20), see [ ]. The divergence of 2'^,fc(?) for ^ oc makes the finding of high 
real mode numbers, especially at large fc, extremely difficult. We stabilize the search by 
using an analytic approximation of the highest real mode number z^max(fc5 o) = niax{z/ : 
Ziy,k{Q) = 0} for q sufficiently far away from 0. This is achieved by using the asymptotic 
expansion of Bessel functions [ ] in the transitional regime yielding 



. 23a2 
ao + ai exp k 



+ 0(fc-3), (23) 



where 



^1 = — = 0.323685 , ao = 2.3381074 . (24) 
2Ai (-ao) 

The constant ao is the negative first zero of the Airy function, Ai(— ao) = 0. The exact 
implicit formula Zy^k{q) = for mode numbers at given k and q ^ has an interesting 
simple first order approximation reading 



/co(n, q) J \ nn 

which is asymptotically exact in two independent limits: u = fixed, ^ oc, and 
k = fixed, oo. The relation (25) represents a useful approximation of mode 

numbers and is to our knowledge a new uniform approximation of modes in a bend. The 
expression /co(n, q) is the n-th zero of Zo^k{Q)^ which can be easily found numerically. In 
the limit of large n, where we can use /co(n, q) ^ nn/a^ the relation (25) in simplified to 

{kaf - {lylogqf ^ {Tinf . (26) 

The validity of this formula is illustrated in figure 2, where we compare mode numbers 
obtained from the approximate relation (26) with the exacts ones. The highest real mode 
^max for small wave-numbers k < u can be approximated using equation (25) as 



TT / / k ^ ^ 



In practical applications it is important that below the wave- number k\ow{q) = ko{l^q) 
there are no real modes. In wide bends with q ^ one can expand the cross-product of 
Bessel functions around q = and obtain the formula 

hoM)=bo + T7^ + 0{q'), (28) 

where bo is the smallest zero of the Bessel function Jo{x)^ Jo{bo) = and 

bo = 2.404825558 , bi = -^^7^ = 1.54288974 . (29) 

In narrow bends, where a(= 1 — g) ^ 0, we can use standard stationary perturbation 
theory [ ], see equation (B.2) in Appendix B, to approximate k\ow{q)' By introducing 
the matrix elements 

1 ^ sinfymx) sinfyrmx) q , 

Vnm = -7^ dx \ ' ^, 7=7^, 30 

2 Jo (x + 7)2 1-q 
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(a) (b) 

Figure 2. Two different representations of mode numbers u at the corresponding 
wavenumbers k in the bend with q = 0.6. In (a) we plot k vs. where we separately 
discuss real and imaginary mode numbers depicted on the right and left side of the 
abscissa, respectively, and in (b) we plot k'^ vs. v'^. The dashed lines represent solutions 
of equation (26). 



a 



we can express the lowest wave-number as 

-I 1/2 

(31) 

1 — L" 

1>1 

The formula (31) has a simple first order expansion in a = 1 — g 

koM = --^ + 0{a') . (32) 

From the expression (32) we learn that the lowest wave-number at which real modes 
exists increases with increasing q and converges to ^. 

2.2. Numerical evaluation of mode functions in a bend 

The mode functions in the bend at a given wave-number k and inner radius q 
are proportional to cross products of Bessel functions Zjy^ki"^) (12), where the order 
parameter u takes values from the set Aik^q- Because of the symmetry Z_^^/e(^) = Z^^ki"^)^ 
we only consider mode numbers from the set Aik,q,+ = {^n ' ^ ^n+i ^ ^ N}, which 
are ordered by decreasing square. To illustrate the basic properties of mode functions, 
we plot in figure 3 the functions Zj^^^i^) for real and first few imaginary mode numbers 
for q = 0.6 and some low wave- number k. We see that the first mode function Zj^^^ki"^) 
has no zeroes on the interval r G (?, 1) and each consecutive mode function has one 
additional zero. In the following, we present formulae and numerical recipes for stable 
calculation of mode functions, where we assume that the mode numbers are given. 

The Bessel functions for real orders are well implemented in the currently available 
numerical libraries, i.e. SLATEC [22]. From the definitions of Bessel functions it is not 
unexpected that we encounter problems at evaluating Zj^^k for large wavenumbers k ^ 1 
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Figure 3. Plots of Zj^^^kij) for first 12 mode- numbers Vn ^ (^n ^ ^n+i) 

calculated at the wave-number k = 78.5398 (TVb = 10) and inner radius q = 0.6. 



and high real orders u 1. In order to overcome these problems we apply the stable 
forward recursions in the order parameter [ ] written as 

= Pi^-i - 2m(^<7m + yr,^) , (33) 

= (2/x(/x + 1)?/^ - l)r^ + 2/x(/x + l)xyq^ - + + l^xp^ , (34) 

r^+i = (2/x(/x + l)x^ - l)q^ + 2/x(/x + l)xyr^ - (m + l)xp^-i + jiyp^ , (35) 
where we write x = {kr)~^ and y = and define the following symbols 

= J,{kr)Y,{k) - Mk)Y,{kr) , (36) 

q, = J,{kr)Y'^{k) - J',{k)Y,{kr) , (37) 
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r, = J',{kr)Y,{k) - J,{k)Y'^{kr) . (38) 

The initial conditions for the recursion, at low orders, are calculated using standard 
routines and the expressions for the derivatives of Bessel functions, e.g. a relation valid 
for any Cylindrical function: = (Cj^-i — Cj^j^i) /v. However, we encountered a problem 
at high wavenumbers and low orders due to the lack of precision in SLATEC routines. 
Therefore in that regime we use the Hankel approximation [19] to evaluate Z^^]^ 

Z,Jr) = [(AC + BD) sm(k(l - r)) + (AD - BC) cos(fc(l - r))] , (39) 

where we write A = Pjj{kr)^ B = Qiy{kr)^ C = P^{k) and D = Q^{k)^ which are 
expressed in terms of the asymptotic series: 

(^- 9) (^-l)(^-9)(/i-25)(^-49) 

2\{%zf ^ 4!(8^)4 ^^"^^ 

and 

which we sum up to /^ax = L^/2j + 1 terms. There is another difficulty occurring at high 
wave-numbers, which can not be corrected. The first few real modes G Aik,q,+ scale 
linearly with the wavenumber k and functions ^^^^/^(r) diverge with increasing order u. 
Consequently, the values Z^^k{(l) are exponentially sensitive on the precision of the first 
few mode numbers v^. 

\Z,,+5u,kM ^ sinhd loggi 5p) . (42) 
This problem can not be solved completely, but only partially corrected by manually 
setting the values of ^^^^/^(r) to zero around the inner radius. This can be done without 
any real loss of precision, because the mode functions are localised near outer radius. In 
practice we calculate the left side of equation (42) in a finite range arithmetic [— m, m] 
with the maximal number m (e.g. m ^ 10^°^ in double precision). By considering that 
together with the known property z^max ^ k (23), we find that in practice our modal 
approach breaks down above some wave- number fcbreak ^ logm/| logg| and consequently 
bounding the number of open modes in our numerical analysis below 

a logm 

^Vo,break ~ " "H T ^ l^dj 

TT I logg| 

where we used the relation A^o ^ ka/n valid for narrow bends and high wave- numbers. 

The numerical evaluation of Z^^^ki"^) at imaginary orders z/ G iM is almost 
unsupported in currently available numerical libraries. Therefore we have developed 
procedures for their evaluation ourselves and give here a summary of our work. For 
each regime of order parameters ly = iy and wavenumbers k we use a different strategy 
to evaluate Zi^^/e(r) in order to achieve an optimal precision control and a CPU time 
consumption. The formula for the cross-products of Bessel functions (12) takes for 
imaginary orders a simple form 
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At small wave- numbers k or more generally for /c <C we use the Taylor expansion of 
the Bessel function [19] and rewrite equation (44) into 

Ziy,k{r) = {r-'yu,yikX^{kr)} , (45) 

where we use the series 

The series (46) is summed up to the index /^lax = + |^^p/4j , where e is the 

desired accuracy of the expression. 

At higher wave-numbers k and orders y ^ k we combine the backwards recursion 
valid for Cylindrical functions 

ai_i{z) = ^^^^^^\ i{z) - a/+i(z) , (47) 

with an appropriate normalisation formula for ai{z) [ ] and thereby obtain an 
expression for u\y{k) (45) given by 



Uiy{z) = ao{z) 



^ {2l + z){l + z,l + z) 

.1=0 



Ci2l{z) 



(48) 



The terms a2i in the series (48) are given with the recursion (47) started at the index 
^max = 2[(x + l)/2j with initial conditions ai^^^ = £ and a/^^^+i = 0, where constant £ 
is the smallest number supported by the CPU architecture and x is determined by the 
equation 

I log^l — 1 + - log(l + y^) — y arctan y = 
\ f z \ y 

-X log(x^ + y'^) — X yiog - + 1 j — X arctan — . (49) 

The later equation (49) is meaningful only if the right side is positive, yielding that 
presented approach with the iteration formula is valid only for orders y below some 
value scaling as 0(|log^|). By increasing the wavenumber k further up and keeping 
orders small y ^ k we can use the Hankel approximation (39) with the order parameter 
u = iy. The asymptotic series (40) and (41) in the expression (39) are summed at least 
up to /max = 1^/"^ + \/^^^^^-yy^\ terms. At large enough wavenumbers k and imaginary 
orders y > k we can make use of the Debye approximation of Bessel functions [23] and 
write 

2 ^> \G{a, i)G\(5, ^)e4^-C-K™h|-arcsinh^)] I 
ZUk, r) = ^ , , ^ ^ , (50) 

with substitutions 

e = vW^, C=Vy' + (kry , a-i = i+ /?-i = i + L_L. (51) 
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The expression y) in formula (50) is given in the form of an asymptotic series 

G,..,) = f;(ziM. (52) 

m=0 ^ 

where polynomials Vjy^it) are generated by the following recursion 

1 ]^ fc + l k 1 

Vk+l{t) = -{l-t){kVk{t) + 2tv'k{t)) + —t-— / {l-bT)T — Vk{T)<\T. (53) 

2 16 Jo 
The first few Vmit) read as 

/\ - /\1 5 3 TT 385 9 

v^it) = 1, ^i(t) = 1, ^2(i) = 1 + r , . . . . (54) 

^ ' ^ 8 24 ' ^ ^ 128 576 3456 ' ^ ^ 

The formulae (39, 45, 48, 50) enable a stable high precision calculation of the mode 

functions in the bend at imaginary mode numbers. 

2.3. The overlap of mode-functions in different geometries 

The main ingredient in the modal description of the scattering are the overlap integrals of 
mode functions in the straight waveguide and in the bend. These overlap integrals "tell" 
about the compatibility of both scattering regions and are discussed in the following. 

The cross products of Bessel functions with order parameter u G Aik,q,+ at given 
wavenumber k and inner radius q form a set of functions 

Zk^q = {Zr.,k{r) : ly e Mk,q,+, T G [g, 1]} , (55) 

which is complete in L2[q^ 1] and orthogonal w.r.t. the weight function w{r) = r~^. The 
later is derived in Appendix A. The orthogonality relation for Zj^^^ ^ 2k^q reads [24] 

k 

2^ 

The line separating the bend and the straight wave-guide will be called the cross-section 
of our open-billiard. On the cross-section we define two different scalar products denoted 
by (•, •) and (•,•), and written as 

(a, b) = f -a{r) b{r) , (a, b) = f dy a{y) b{y) . (57) 

Jq ^ Jo 

Let us now introduce modes at some fixed wave-number k and inner radius q for 
different regions of the open billiard. In the bend, mode numbers Up and normalised 
mode functions Up{r) read as 

u,eMk,,,+ , U,{r) = ^=^d^= , PEN, (58) 

where we order the mode numbers so that v"^ > Up_^_^ , and in the straight leads connected 
to the bend we have mode numbers Qn and mode functions Un{x) defined by 

- , Un{y) = sin (^ni/) , n e N . (59) 



^i.+i,/c(l)^^(l) - qZ^+i^k{q)^^^{q) 



(56) 
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The modes with real and imaginary mode numbers are caUed open and closed modes, 
respectively. The number of open modes in some geometry is denoted by A^o- The overlap 
integrals of mode functions are given by 

Anp = {Un, Up) , Bnp = {Un, Up) , (60) 

where we use the relation r = q + y between the coordinates. In figure 4 we show a 
density plot of the matrix elements, in log scale, namely log \Anp\ and log at two 

values of inner radii with the same number of open modes A^o in both geometries. We 



g = 0.2 q = 0.9 




1 50 1 00 150 300 250 3O0 1 50 100 150 300 250 300 




-10 -9 -8 -7 -6 -5 -4 -3 -2 -1 



Figure 4. Density plots of the matrix elements log;LO \^np\ (top-row) and log;LO l^npl 
(bottom-row) for q = 0.2 (left) and q = 0.9 (right) as indicated in the figure. The 
number of open modes is A^o = 100 and the number of all considered modes is = 300. 

see that the matrices A^p and B^p have a similar form for all q and k. This is starting 
at small indices with a wide area of high values of matrix elements that squeezes to 
almost a single intensified point at A^o again spreading in a triangular shape with 

increasing indices. The parameter q has a strong infiuence on the shape of the area with 
high intensities in A and B. In the case of small values of q in contrast to larger the area 
of high values in matrices A and B covers almost the whole open-open block of indices 
and with crossing of the narrowing at n,]9 A^o spreads faster with increasing indices. 
The shape of matrices A and B is similar therefore in the following we only show results 
for the matrix A^p. We found numerically that the area of high intensities in matrices A 
and B scales with the number of open modes A^o as Anp^ B^p ^ F{n/No^p/No)^ where F 
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is some a well behaved function. We demonstrate this by plotting the matrix elements 
\Anp\ in relative indices n/No and p/N^ for different numbers of open modes A^o shown in 
figure 5. In addition, we find numerical evidence that intensities of matrix elements A^p 



No = 100 No = 200 No = 400 




-10 -9 -8 -7 -6 -5 -4 -3 -2 -1 



Figure 5. Density plot of matrix elements log^g \^np\ at inner radius q = 0.6 and 
different number of open modes A^o = 100, 200, and 400 as indicated in the figure. 



and Bnp in the region of open modes and in the region of closed modes scale differently 
with A^o 

1 



open modes : |A 
closed modes : \A, 



np\ -) 



\B. 



np \ -) 



np I 



~ ]\[2 



n p 

iVo'iVo 



n p 



(61) 
(62) 



where F^ax and Ftaii are some well behaved functions. Taking into account these 

phenomenological findings enables a better precision control of scattering calculations. 

1 

In figure 6 we see that F^ax is an envelope function for maximal values of A^o^ \Anp\ and 
that Fmin can be chosen to fit the tails of A^^l^^^l. We prove the scaling relation (62) 




Figure 6. The cuts in the index space of matrix A^p at fixed p/N^ = 0.5 and different 
number of open modes Nq = 50, 100, 200, and 400 at inner radius q = 0.6. 
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by using an asymptotic approximation of the mode function in the bend 

Tip 



sm 



logr 



loggi V^^SQ' 

which yields the foUowing asymptotic behaviour of matrix elements 
6 |logg|i 



A 



np 



71' 



{q^-iy + (-1)"+') np-^ + 0{p- 



A 



np 



(l-<?)^ 



n 



fixed, p ^ oo 



(63) 



(64) 



2 I 1 I 3 

^ |l0gg|2 



p fixed ^1, n ^ oo . (65) 



The off-diagonal diagonal elements A^p and B^p decay algebraically with increasing 
index. The pre-factor of the decay is decreasing with increasing q and is singular at 
g = 0. This means that for large enough q it is possible to approximately express the 
open modes of the bend solely in terms of open modes of the straight waveguide and 
vice versa, as indicated by the relation (61). 

From the definition of matrix elements A^p and B^p (60), and completeness of the 
mode functions at given k and g, it follows that A and B are transition matrices between 
the sets of mode functions in different regions, 

^n{x) = ^ Bnp Up{r) , Up{r) = ^ A^p Un{x) , (66) 

yielding the relation 

AB^ = A^B = id . (67) 

In practice we work with finite sets of modes, where the identity (67) cannot hold 
exactly. In figure 7 we plot AB^ for different inner radii q and fixed A^o = 100. The 



g = 0.2 



q = 0.6 



■ 
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Figure 7. Density plot of matrix elements log^g \ {AB^)nm — ^nm\ for q = 0.2, 0.6 
and 0.9 (from left to right). On the abscissa and ordinate we plot indices n and m, 
respectively. The number of open modes is A^o = 100 and the total number of modes 
is 300. 



mismatch from the identity (67) on some sub-set of indices starting at the origin n = 1, 
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or p = 1, increases with decreasing inner radius q. This means that the numerical 
calculation of scattering for smaller q should be less accurate at finite dimensions. The 
discrepancy between AB^ and the identity on (truncated) finite dimensional spaces is 
strongly non-uniform in indices. Before going into practical aspects of this problem, we 
examine the convergence of matrix elements {AB^)nm to Snm with increasing number 
of all considered modes N = Nq + Nc Sit fixed A^o, where A^c is the number of closed 
modes. An example of such convergence is shown in figure 8. The convergence proceeds 



A^c = 1 A^c = 100 A^c = 500 




1 io pn ^ ^n rn Til fin sn ' m 1 pn 1:1 Rn Rf^ ' m ■ m n:i ^ fa ' fln Tvn 1 IiTTi ?m vni nnn r£Y\ ?ff] 
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Figure 8. Density plot of matrix elements log^g \ {AB^)nm — ^nm \ for different numbers 
of closed modes TVc = 1, 100 and 500, as indicated in the figure, at inner radius q = 0.2 
and A^o = 100 open modes. The labels on the abscissa and the ordinate are matrix 
indices n and m, respectively. 



by the standard scenario, where the agreement between [AB^]np and 5np propagates 
block- wise from low to higher indices with increasing A^. The propagation is slow due 
to the triangular shaped area of high intensities in the closed-closed modes block of A 
and B. The speed of propagation of accuracy to higher indices increases with increasing 
inner radius q. 

The SVD decompositions [ ] of matrices A and B is useful for improving and 
stabilising the scattering calculations and will be used in the next section. From 
definitions of transition matrices (60) and completeness of mode functions we obtain 

{AA^)nn' = {un, r Un') , (^^^W = (^p. ^ . (68) 
{BB^)^^> = {un, r-' Un') , {B^B)pp. = {Up, r'^ U^^) . (69) 

which we use to bound the image 

for ||a||2 7^0. (70) 

These results together with AB^ = id can be used to determine the form of the SVD 
decomposition 

A = UEV^ , B = f/E- , E = diag {a, e [^q, l]},eN ' (71) 



Hallo 



|-D^I|2 

Hallo 
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where U and V are orthogonal matrices. We show here an example of the SVD 
decomposition of finite dimensional matrices A and 5 at g = 0.2 and A^o = 100. In figure 
9 we show singular values and in figure 10 we show density plots of the corresponding 
matrices U and V ^ where the inner indices are ordered by decreasing magnitude of 
singular values. We see that the relative dimension of the space, which violates 
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i/N 




0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 
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Figure 9. Singular values of matrices A (left) and B (right) for different number of 
closed modes = 100, 200, 300, 400 and 500 at q =0.2 and = 100. 





Figure 10. Density plots of matrix elements \Unm\ (left) and | Kirn I (right) of the SVD 
decomposition of matrices A (top row) and B (bottom row) at inner radius g = 0.2 
with No = 100 open and Nc = 100 closed modes. Labels on the abscissa and the 
ordinate are indices n and m, respectively. 
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the bounds of the singular spectra (70), converges with increasing space dimension 
N = No + Nc^ where A^c is the number of imaginary modes. In the presented case q = 0.2 
the relative dimension is around 20%. It is important to see that vectors in V and U 
corresponding to singular values, which violate the bounds, have non-zero components 
only at closed modes. We conclude that due to a finite dimensional representation of 
matrices A and B we have deviations from the infinitely dimensional case only at high 
laying decaying modes that span a space of almost fixed relative dimension for some 
value of q. 



3. The scattering across a bend 

In this section we solve the on-shell scattering of a non-relativistic particle across our 
open billiard using the modal approach, initiated in the introduction. The scattering 
is discussed at fixed wave-number k and inner radius q. The control of the precision 
of scattering calculations is studied in detail. In the second part of this section we 
investigate some interesting physical scattering properties of the bend. 



3.1. The scattering matrix of a bend 



The scattering matrix S [ ] is a linear mapping between the incoming and outgoing 
"waves" with respect to our scatterer. We reorganise the expansion coefficients 
of the wave function over different regions (9), (10), (11) into the incoming 
contributions denoted as = (a^, a^, . . . , 6^, 63^, . . .) and outgoing contributions 
'^out = (^^r^ ^2 ^ • • • ^ • • •)• Then the scattering matrix S can be defined as 

Svi^ = Vout . (72) 
The S-matrix in our case has a simple symmetric block form 
R T 



S 



T R 



(73) 



with R and T being the refiection and the transmission matrix, respectively. By 
reordering of rows and columns in S so that the matrix elements concerning open 
(subscript o) and closed (subscript c) modes are separated and grouped together we 
obtain the matrix 5, reading 



^^^^ S S 

We find that blocks of matrix S (74) obey the generalised unitarity 
following relations 

id , 



^00 ^lo 

i 5co5]o 



^00^00 

^OC 5 
Sen 1 



i ^oc^oc 



(74) 
defined by the 

(75) 
(76) 
(77) 
(78) 
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The relations (75 - 78) result from the probability current conservation, which is also 
equivalent to the condition AB^ = id. Due to the time-reversal symmetry of the physical 
problem the scattering matrices S and S are symmetric 

= S , = S . (79) 

The block symmetry of the scattering matrix (73) simplifies its calculation. We may 
consider individual incoming waves represented by the following wave function ansatz 

^(r) = e+(r) + J] e" (r) , r e I^a , (80) 

m 

V'(r) = E ^tn + ' r G , (81) 

P 

^(r) = J]e+(r)r^„. r e fie , (82) 

m 

The continuity of ip{r) and its normal derivative on the connecting cross-sections between 
regions f^A,B,c determines the matrix elements Rmm T^n and A^^ and yields the following 
system of matrix equations 

id + i? = M(A+ + A") , T = M(j^A+ + T'^A') , (83) 
id - i? = A^(A+ - A-) , r = A^(j^A+ - T-^A-) , (84) 
M = G^AV-^, N = G-^BV^, (85) 

which we write using diagonal matrices V = diag{z/^}^^i^, G = didig{gn}neN and 
T = exp(il/), and transition matrices A and B (60). The elimination of matrices A^ 
from equations (83) and (84) yields the blocks of the scattering matrix reading 

T={C- TDC-^TD)-^T{C - DC-^D) , (86) 

R={C- TDC-^TD)-\TDC-^TC - D) , (87) 

that we express by using the following auxiliary matrices 

C = M-^ + A^-^ = V^B^G-^ + V-^A^G^ , (88) 

D = M-^ - N-^ = V^B^G-^ - y-^A^G^ , (89) 

where we take into account the relation AB^ = id. The presented form of the matrix T 
(86) and R (87) is chosen in order to increase its numerical stability i.e. minimising the 
use of inverses and avoiding direct computation of 



3.2. Numerically stable scheme for scattering matrix calculation 

A bend on a straight waveguide is a paradigmatic example for testing numerical schemes 
and ideas on how to accurately calculate the scattering matrix. In particular, the high 
curvature case q ^ turns to be highly non-trivial. Here we give a simple and stable 
procedure to obtain the scattering matrix with a clear precision control for practically 
all curvatures. 

The scattering across a bend of angle (3 and inner radius q back to asymptotic 
region at some wavenumber k is described by the scattering matrix S{(3) (73), which is 
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composed of the reflection matrix (87) and the transmission matrix r(/?) (86). In 

practice, we work with flnite dimensional matrix approximations, denoted by 

where > A^o is the number of modes used in the asymptotic regions. The main 
objective is to construct these flnite dimensional matrices Rn{P) and Tn{I3) so that: 

(i) calculations are numerically stable and precise, 

(ii) Sn{(3) satisfles the time reversal symmetry (79) and the generalised unitarity 
relation (75 - 77), 

(iii) the sub-block of i?Ar(/3) and T/v(/?) of dimension A^' G (A^o,^] is calculated with 
controllable accuracy, where A^o = L^^/^J is the number of open modes in the 
asymptotic region. 

The recipe to achieve these assumptions may be separated into two parts. In the flrst 
part we cure the numerical instability caused by the maximal element r(/?) = exp(/?z/7v) 
in ^7v(/?), which are exponentially diverging with increasing A^. This is achieved 
by separating the bend into 2^° identical subsections of angle (3' = 2~^^f3 so that 
^(/^O^num ^ 1, where £:num is the numerical precision e.g. ^num = 2~^^ in double precision 
floating point arithmetic. The number no is calculated as 

log 6:nun 



no max < 0, 1 + 



log2 



(91) 



The scattering matrix S'Ar(/3') of a small subsection of the bend can be calculated in 
a very stable way and up to a high precision. By concatenating scattering matrices of 
subsections together with the recursion 

5^(2--+^/?) = 5^(2--/?) 5^(2--/?) . (92) 

we obtain the scattering matrix of the whole bend Sn{P). The symbol denotes the 
operation for concatenating scattering matrices associated to scatterers on the waveguide 
and is deflned in Appendix C. 

In the second part we are discussing the problem that SN{f3^) diverges with 
increasing A^, which is due to violation of the identity AB^ = id (67) for flnite truncated 
transition matrices and B^. We eliminate the problem by deforming A^ and B^ 
so that they are non-singular and exactly fulfll the condition A^Bj^ = id. We make 
the SVD decomposition of the truncated matrix A^ = UnT>nV^^ modify its singular 
values Sat = diagla^}^;^ to Sat = diagja^l^^^ so that they flt in the bounds obtained 
for inflnitely dimensional case (71) 

a.= h' ■ '^^^^[^^'l] , (93) 
I (J* : otherwise 

and again generate both matrices 

An = Un^nV^ , Bn = C/ivS^V^ . (94) 
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The same procedure can also be done using SVD decomposition of the matrix slS 
a base for generation of both deformed matrices and B^. The value of a* > can 
be chosen arbitrarily, but the most elegant choice is a* = 1. By using matrices 
and Bn (94) instead of A^ and B^ in SN{f3^) and consequently in SN{f3) these become 
generalised unitary with a well behaved and physically precise limit ^ oc at least 
on the sub-space of dimension N\ We first check this by discussing the precision of 
transition between the modes in the bend and in the straight waveguide on the sub- 
space of dimension < N. The error of transition from the asymptotic region (infinite 
waveguide) into the bend is quantified by 

es^iy{N^ N^) max |(u 

and for transition in the opposite direction by 
£b^s(iV, iV')=max|([/j 

with L = {1, 2, . . . A^'}. The introduced transition errors (95) and (96) measure the 
violation of the identity AB^ = id on the subspace A^^, when working with the finite 
number of modes A^. It is expected and supported by our numerical studies that the 
errors vanish in the limit A^ ^ oc at fixed A^' and other parameters. In figure 11 we show 
transition errors as a function of A" at a fixed A^^ for two values of A^o- The errors decrease 



N 



N 



,)| = max 



(95) 



N 



N 



^ \Un){Un\ - id 



\Ur] 



n=l 



max 

p,rEL 



n=l 



nr-^np 



pr I 



(96) 




Figure 11. The transition errors Ss^hi^o + Nc, N') (top) and £:b^s(A^o + A^c,A^O 
(bottom) calculated on the sub-space of dimension N' = Nq -\- 10 vs. the number of 
closed modes for = 10, 100 (left, right). 



with increasing A^c down to a certain plateau, which is determined by the precision of 
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mode functions. The most problematic in the precision are the highest few open modes 
in the bend. A conservative estimate for the plateau is around 10~^. The transition 
errors (95) and (96) increase with decreasing q indicating that we need more modes to 
achieve equally small error as for higher q. We did not find any analytic approximation 
for transition errors. Therefore we numerically estimate the minimal dimension of the 
functional space A^t = + Aq needed for transition errors to be smaller than some 
defined as 

At(A^ e) = min{A : 5s-b(A, A^ < s and 5b-s(A, A^ < s} , (97) 

where A^ = N^+Nq is the dimension of the observed sub-space. In our numerical analysis 
we set £ ^ 10~^. We check the convergence of Rn{P) and Tn{(3) with increasing A on 
some fixed sub-space of dimension N' < N . The matrices Rn and T/v are calculated using 
the method of dividing the bend into sub-sections together with deforming the transition 
matrices. The convergence is measured through the relative difference of matrices at 
subsequent changes of the dimension A 

^ / ^T,^ \\Rn -rl{[^) - Rn{(^)\\n' ( . wx _ ||7V+i(/3) - Tn{I3)\\n' .q^x 

'■''''•'^> = — WimU' — ' — ' 

where we introduce a matrix norm ||^||m = niax^j^[i^M] on the sub-space of 
dimension M. The expressions (98) give upper bounds for the deviations of matrices 
from their asymptotic forms 

CXD 

\\RooW)-RNm\N'<Ci{N') en{M,N'), (99) 

M=N+1 
OO 

\\Too{P)-Tr,{mN'<C2{N') ^t{M,N'), (100) 

M=N+1 

with expressions Ci{N') = maxAr>Ar/ ||i?Ar(/?)||Ar^ and C2{N') = maxAr>Ar/ ||T/v(/3)||Ar^ 
which are of the order of magnitude 1. We have numerically studied quantities 
^r,t(A, A^) as functions of A at fixed A^ and the results are shown in figure 12. In 
the case A = A^ = Aq we are talking about the open-open block of the scattering 
matrix, which is sometimes called the semi-quantal approximation [ ]. We see that £:r 
and decrease with increasing total number of modes A = Ao + Ac in a similar fashion 
as transition errors down to some plateau around 10~^^. The plateau is almost equal to 
the machine precision, which is surprisingly better than the transition errors. 

The presented method for the calculation of the scattering matrix and its accuracy 
control works well for all g > 0. Nevertheless, a treatment of high curvature cases g ^ 
are difficult as we need to consider a large number of closed modes to reach a sufficient 
precision of the scattering matrix. However, this is feasible to achieve by our method in 
a stable and controlled way. 

3.3. Quantum transport across the bend 



The scattering matrix of the open billiard describes the stationary quantum transport 
of a particle over the bend. We discuss here the most important and obvious measures 
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0.01 




Figure 12. The measures of convergence 6:r,t(^o + N^^N^ + 10) (top, bottom) as 
functions of A^c for = 10, 100 (left, right) at various inner radii q as indicated in 
the figure. The bend is of angle /3 = tt. 



of the transport, which are the reflection probability and the Wigner- Smith delay time. 
Since the classical particle cannot scatter back [ ] we are particularly interested in 
the reflection probability as a genuine quantum (wave) property. This is a measure of 
quantum tunnelling between the two classically invariant components of the phase space 
corresponding to right and left going waves. The scattering of an incoming Gaussian ray 
over the bend is illustrated in figure 13. The ray follows the classical trajectories, but 
as it is of finite width, its parts are scattered differently when hitting the curved wall. 
The parts of the ray travel different lengths and interfere among themselves. 

We are discussing the scattering over the bend at some fixed wavenumber k and 
inner radius q. The scattering properties are contained in the transmission matrix T 
(86) and the refiection matrix R (87). The wavefunction over the asymptotic region is 
described in modes, where N > Nq = [/ca/vrj . We consider an incoming wave coming 
to the bend from the left side written in the asymptotic region as 

N 

^i„(r) = J]a„e+(r). (101) 

n=l 

By introducing a vector of complex coefficients a = {an}n=i we can write the transmitted 
and refiected probability fiux, and j^r respectively, in an elegant form 

Jr = a^Ua , Jt = a^'Sa , jo = a^a = Jr + Jt , (102) 

where we introduce matrices IT and S calculated from open-open mode blocks of R and 
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Figure 13. Scattering of an incoming ray of Gaussian shape with unit probabihty flux 
at wavenumber /c, which supports 100 open modes, and for two inner radii q = 0.2, 0.6 
(left, right). 



n = RlRoo , S = r^t^Too , where H + S = id . (103) 

The average transport properties are given by the first and the second moment of the 
probability currents jr^t averaged over an ensemble of incoming states a. The ensemble 
represents states (vectors a) uniformly distributed over the 2N-dimensional sphere of 
radius jo [18, 29]. The average probability currents are given by 

n = ^^^ = ^tT{u}, r=l-7^, (io4) 

and the standard deviations of probability currents (giving fluctuations within an 
ensemble) are written as 



2 _ (Or - Or)^)^)^ _ 1 

(Jo — 



tr{n^}_^, 



In the following we thus discuss the average reflection TZ and the dispersion of reflection 
a^. An approximation of the transmission matrix T can be determined from the semi- 
classical calculations, whereas for the reflection matrix R it can not, as the reflection in 
the bend is a purely quantum phenomenon. The gross structure of matrices R and 
T, similarly as of matrices A and 5, does not change signiflcantly with increasing 
wavenumber k. In flgure 14 we show the density plot of matrices R and T with A^o = 100 
open modes. The high probabilities in the matrix T have a classical correspondence, 
which is revealed through the calculation of the classical scattering matrix Td^ss [ ]• 
Both, the classical and the quantum transmission matrices feature similar patterns, 
but due to the quantum interference, we can not establish a clear correspondence. In 
the matrix T we have a large area of high values so we can expect that transmission 
probability of individual modes should be high. It is important to notice the area in the 
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reflection matrix of high intensity is concentrated around the last open mode with the 
index A^o- 

Rn Tji Tclass 




1 £0 40 6 80 1 00 120 140 1 £0 40 60 80 100 120 140 1 10 20 30 40 50 60 70 80 90 100 
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Figure 14. The density plot of the scattering matrices log^o |(^Ar)nm| and 
log;LO \ {TN)nm\ and the classical analog of the later log^g |(^ciass)nm| calculated at inner 
radius = 0.6, wavenumber k = 100.5^ and with the total number of modes A/' = 150. 



At the so called resonant wavenumbers fc^ = (m G N) a new open mode appears 
in the asymptotic region and causes a strong increase in the reflection matrix elements 
at open modes with high indices. This is demonstrated in flgure 15, where we show the 
scattering matrices R and T around the highest open mode calculated calculated at 
k ^ fcioo and at g = 0.6. The changes are centred around the index A^o and signiflcantly 
influence the average transport. 

In order to clarify the contributions to the total reflection we plot in flgure 16 
the reflection probability of individual modes H^^ for wavenumbers near and far from 
the resonance. We see that the highest open mode has the strongest reflection and 
the reflection probability of "all" modes increases at resonant wavenumbers k = k^^. 
In particular the reflection of highest open mode is almost perfect nAro,Aro 1. In 
the vicinity of the resonance we could effectively approximate the average reflection 
diS TZ ^ Un^^n^/Nq (see flgure 17). The resonant wavenumbers k^^ are important 
markers for anomalously strong reflection. This is illustrated in flgure 17, where we 
plot the average reflection 7^ as a function of the wave- number around k = fcioo- We 
see the average reflection TZ has a strong sharp maximum at resonant wavenumbers and 
decreases in an irregular oscillating manner with increasing wavenumber until crossing 
the next resonant wavenumber. The frequency of irregular oscillations increases with 
increasing wavenumber. From numerical results we see that TZ decreases with increasing 
q. In narrow channels a ^ at large wave- numbers we showed with a perturbative 
approach (see Appendix D) that 

^ - ^ , (106) 



The bends on a quantum waveguide and cross-products of Bessel functions 25 



before middle after 




-9 -8 -7 -6 -5 -4 -3 -2 -1 



Figure 15. The density plots of the reflection matrix log^g |(^Ar)nm| and the 
transmission matrix log^^g l(^A^)nm| around the point (n, m) = {Nq^Nq) at 
wavenumbers k = (100 — 10~^)^ (before), k = 100;^ (in the middle) and k = 
(100 + 10~^)^ (after the resonance). The solid lines denote indices n = Nq and m = Nq. 
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Figure 16. The diagonal matrix elements Ilnn at wave- numbers k = as indicated 
in the figure, and fixed inner radius q = 0.6. Note that A^o = M • 



which is confirmed numerically. The resonant behaviour around the resonant 
wavenumber can be partially explained by neglecting all open modes expect the one 
with the highest index A^o- Such system can be treated as an independent Id scatterer 
{d = 1) with the reflection and transmission matrix elements reading 

_ sm(/?^7vj sin(i/x^J 
sm(/^z/7v, + i/xatJ sm(/^z/7Vo + i^iVo) 
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with the phase shift 

IJ^^ = 2 arctanh (K^J , = G [0, 1] . (108) 

At = /caTo, the mode number g^^ ^ fi^^ and the phase shift fi^^ become zero yielding 
a perfect reflection in a Id scattering model, i?id = — 1 and Tid = 0. This treatment 
is meaningful, because the matrices A and B are approximately diagonal at (Nq^Nq) 
with an algebraic decay of matrix elements when we move away from the diagonal. If 
the modes were strictly independent we would have n7Vo,Aro = H^idlP, but the algebraic 
tails in matrices A and B make this solution to hold only as a rough approximation as 
can be seen in flgure 17. In flgure 18 we study measures of reflection TZ and ctr over a 




k( 1 -q)/7l - 1 00 1 0O(r-rO), r = k(1 -q)/% 



(a) (b) 

Figure 17. Reflection a bend of an angle ^5 = tt as a function of the wavenumber k 
in the region of a reflection resonance k = 101^. In (a) we show the average reflection 
measure for different values of q (indicated in the figure), while in (b) we show a zoom- 
in around the resonance and compare the average reflection and its standard deviation 
with the reflection in the last open mode (inset shows even a more very narrow region 
of the resonance). 



larger range of wave-numbers k. From flgure 18. a we see that TZ strongly oscillates with 
peaks at resonant wave- number kn and its upper bound decreases proportionally to 
as predicted. The numerical results in flgure 18.b indicate that ctr < 7^ and 7^ ^ ctr as 
k goes to inflnity. To get rid of oscillations and get an overall average behaviour of TZ 
and ctr we calculate their cumulative integrals with respect to the wave- number. The 
results are shown in flgure 19 and yield the following dependence 

rk nk 

I UdKoc / a^dK = 0{\og{k)) . (109) 

J ko J ko 

This indicates together with previous conclusions that the reflection measures, averaged 
over small wave- number ranges, indeed scale as 

7^ - aR = 0{k-^) , k^oo. (110) 

It seems that this relation (110) is valid for an arbitrary inner radius q and represents a 
new and very useful information for the study of wave-guides and general billiards that 
include bends. 
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Figure 18. The average reflection IZ and the relative deviation of reflection ctr over 
a larger interval of wave- numbers k at q = 0.6. 




Figure 19. The integral of the average reflection IZ and deviation ctr over larger 
interval of wave- numbers /c of a bend with a radius g = 0.6 and an angle (3 = tt. The 
inserted lines sue fi{k) = 0.01103+0.01022 logA; and /2(A:) = -0.01553+0.006924 log A:. 



Another insight into the scattering properties gives the Wigner-Smith delay time 
Tws [31, 32], which is the quantum analog of the geometric length travelled by a wave. 
In the semi-classical limit, where we could apply geometric optics, r^s is equal to the 
average geometric length of classical trajectories over the open billiard. By using the 
hermitian variant of the lifetime matrix 



°° dik ' 

the Wigner-Smith delay time is defined as 



1 



tr 



°° diA: 



:ill) 



(112) 



where we have used block symmetries of our matrix S (73). can be thought of as an 
average delay time corresponding to particular modes, which are defined as 

1 



m=l 



(113) 
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with derivative defined as {•)' = d/d/c. Numerical results shown in figure 20 point to a 
similar dependence of Wigner-Smith delay time r^s on the wave- number k as the average 
reflection 7^, just the oscillations are smoothed out. The time r^s strongly increases near 
the resonance wave- number k = k^^ due to intense changes in the scattering matrices 
R and T in the area around the index of the newly open mode. Qualitatively we can 
explain the singular behaviour by treating the highest open mode in the resonance 
regions within Id scattering model in which the delay time is given by 

^id _orT*T^ ^p* , _ /?sinh(2/x^Jzy^^ - sin(2/?^^J/x^^^ 

- ^jiidAd + ^id^idl - cosh(2/x^J - cos(2/3z.^J ' ^^^^^ 

The first term in the numerator of equation (114) corresponds to the transmission and 
the second term to the refiection. By slowly increasing the wavenumber across the region 
of the resonance we can notice three different regimes: before, in the vicinity and after 
the resonant wave- number. Slightly before the resonance k < k^^^ a new real mode 
appears in the bend (see formula (B.5)) making the propagation across the bend very 
slow. From the formula (114) we learn that this results in a large transmission time and 
consequently in a large time delay r^g . In the instance of crossing the refiection resonance 
a new mode appears in the straight wave-guide, which causes a square-root singularity 
r^g ^ [k — kN^)~'2 for k > k^^ and its sign is determined by 2(3^^^. This refiection term 
has a short-scale infiuence to the behaviour of the time delay and can enhance or reduce 
its size. Obviously, this is a very non-classical situation. By going further away from 
the resonance wave-numbers the refiection contribution to the time delay is levelled 
by an increasing transmission term due to a very slow propagation of the mode in 
the asymptotic region, which again increases the transition time. So we can experience 
one or two peaks of the time delay in the vicinity of the refiection resonance. Away 
from the refiection resonance the time delay drops even below the classical time. The 
latter we assume is due to refiection phenomena which reduces the classically expected 
phase shift. The presented Id scattering model has only an instructive purpose and does 
not represent any useful quantitative approximation, similarly as was the case in the 
discussion of the refiection. 



4. Conclusions 



We present mathematical, numerical and physical background of the non-relativistic 
2D scattering of a quantum particle on a circular bend connected to infinite straight 
waveguides. We discuss mathematical properties and derive numerical recipes for 
accurate and reliable calculation of the mode functions and the corresponding mode 
numbers in a bend. We take a special care of closed (evanescent) modes in the bend. 
The obtained modal structure and its properties are incorporated in a robust and stable 
numerical scheme for computing the scattering matrix with a controllable precision. Our 
numerical apparatus is applied to the study of transport properties. We focus mainly on 
the refiection, which is a purely quantum (wave) phenomenon since the back-refiection 
of classical rays is not possible. Our study is particularly focused on the possibility to 
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Figure 20. The Wigner-Smith delay time Tws for the bend of angle /3 = tt around 
the reflection resonance at r = 11, 101 (left, right). The solid horizontal lines with 
the labels Tq. 2, 0.6,0.9 represent the classical delay times at g = 0.2,0.6,0.9, which are 
approximately given with the formula Tq = 2.45863 + 2.48696^^. 



investigate the (semi-classical) regime of very large wave-numbers. Some of the obtained 
physical properties of the scattering problem can be explained analytically. In addition, 
we present results on the Wigner-Smith delay in the bend. The obtained transport 
properties can be useful in discussing and predicting properties of open billiards (or 
wave guides) composed of arbitrary combination of bends and straight segments. 
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Appendix A. The symmetry of mode numbers 

We prove the symmetry (15), by changing the variable r = and transforming the 
Bessel equation (3) and the corresponding boundary condition into the equation 



dx^ 



+ {k'e-'^ -iy')Z = 0, Z\,=o^o^, = , (A.l) 
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which can be interpreted as a one- dimensional quantum mechanical eigenvalue problem, 
with the Hamiltonian H and potential V'\ 

a2 ( _h.2 -2x 

V'ix) = \ 
dx^ oo 



HZ = -v^Z , H = 



X e [0, -logg] 
elsewhere 



(A.2) 



Because the Hamiltonian is a Hermitian operator, the eigenvalues are real —u^ G M 
yielding G M U iM. From the form of the potential V\x)^ depicted in figure Al, we 
see that there are only finite number of real mode numbers > and infinite number 
of imaginary mode numbers u'^ < 0. The independent solution of equation (A.2) are 
orthogonal with respect to the measure dx = r~^ dr and so the weight function between 
mode functions Up{r) in the bend is w{r) = r~^. 



oo 



oo 



V'(x) ' 
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-log (q) 
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Figure Al. The analog of the quantum potential in the eigenvalue equation for the 
mode functions in the bend. 



Appendix B. The number of modes in the straight and the bent waveguide 

We discuss the number of open modes N\^{k^q) = card 3?{A^/e,g,+} in the bend and its 
deviation from the number of modes in the straight waveguide = [/ca/vrj. The mode 
numbers in A4k,q,+ continuously shde with increasing k and fixed q from the imaginary 
to the real axis by crossing the point u = 0. This dynamics is depicted in figure 2. This 
means that Niy{k^ q) is equal to the number of zeros x of Z^^xio) up to the value k 

N^{k, q) = card ^{Mk,q,+} = card{x < k : Zo,^(g) = 0} . (B.l) 

By using substitutions U = r~^Lp and r = q + ax the mode problem in the bend (5) is 
in the case u = transformed into a Id stationary Schrodinger equation (7 = ^/(1 — ?)) 

-^ + Vix)^ = e^, V{x) = -[2{x + j)]-\ xe[0,l], (B.2) 

with the eigen-energy dentoed by e = (fca)^. The discrete set of eigen-energies is ordered 
as Cn+i > e^, n G N. By setting V = in expression (B.2) we obtain the mode problem 
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for appropriately rescaled straight waveguide. In the eigen-energies in this case se (t™)^, 
n G N. Then taking into account the (empirical) fact < (vrn)^ < e^+i we can conclude 

< iVb(fc, q) - N,{k, q)<l, Wq e (0, 1) . (B.3) 

This means that at certain k and q we can have in the bend one open mode more, but 
not less than in the straight waveguide. In the semi-classical limit k ^ oo the eigen- 
energies Cn can be obtained using the Debye approximation valid for qk ^ 1. In this 
way we get a relation between the eigenvalues e and its counting number A^b 

27rA^b = (4fc^ + 1)^ - (4(gfc)2 + 1)^ 

- arctan (^(4fc^ + 1)"^^ + arctan (^{4:{qkY + 1)"^) , (B.4) 

which yields with asymptotic expansion in k the expression 

a a(l + q) ^ / / , 



(B.5) 



We see that A^b and are close to each other for high wave-number and and not too 
small inner radius q. 

Appendix C. The method of concatenating scattering matrices 



Here we outline a method to concatenate the scattering matrices [ ] associated to 
scatterers on sectioned wave guides. Let us assume to have two scatterers labelled by A 
and B and with scattering matrices Sx and S'b, respectively. 



S 



A,B 



A,B 



R 

A,B 
R 



(C.l) 



^ (^2Nx2N 

^A,B ' A,B ' 

By combining both scatterers A and B in the order AB we build a "larger" scatterer with 
the scattering matrix S. The matrix S is calculated from matrices 5'a,b by a nonlinear 



operation : C^^^^^^ x 



S = SaQSj, 



r 
t^ 



2^^2iv defined as 



t^ 

R 



r 



(C.2) 



which explicitly reads 

r"- = r\+t^.rU-'t 



A 5 



^B^ ^A 5 



(C.3) 

A 5 t — L^U Ojs^ , (C4) 

where we define L = 1 — r^r^ and V = 1 — r^r^. Note that a bend on a straight 
waveguide can be treated as a scatterer. By combining bends of angles 7^ and 6^ with 
scattering matrices 5^(7^) and S{5')^ respectively, we get a bend of angle 7' + 5' with the 
scattering matrix 5^(7' + 6^). The latter matrix can be obtained from matrices 5^(7^) and 
S{6^) by the formula 

S{i + (50 = S{i) S{5') = S{5') S{i) . (C.5) 
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Appendix D. Perturbative calculation of the scattering matrix for narrow 
bent wave-guide 

We present a semi-classical approximation, for ^ 1, of a scattering matrix 
corresponding to a single bend on a straight wave-guide of width a, as one shown in 
figure 1. Here we are discussing only narrow channels a <C 1, where the influence of 
closed modes on the scattering diminishes. Therefore closed modes are neglected in our 
calculations. We are working at wavenumbers fc, where in all regions of the open-billiard 
the number of open modes is equal. This enables us to write the reflection and the 
transmission matrix in the following simpler form 

R = - iS+]-\C_ + i5_) , T = 2[C+ - iS^]-^ , (D.l) 

where we use the diagonal matrices C = S = ^{^} and G = diagl^'^}^^;^ to 

express the introduced matrices 

C± = G^ACB^G-^ ± G-^BCA^G^ , (D.2) 
S± = G^ASA^G^ ± G-^BSB^G-^ . (D.3) 

We proceed by rescaling the variables to dimensionless form by the following 
substitutions 

y = ^ r = q + ^ tv = ak ^ hn = ag^ \/ ^ — (vm)^ , aVp , (D.4) 

with a new transverse coordinate ^ G [0, 1], and geometric properties being described 
by the parameter a = a/q <^ 1. The transition matrices are then expressed as 

Anp Q^Qnp -) Qnp / = ~~ 5 (-^•^) 

V/oVp(o^(i+«o-^ 

V/oVp(0^(l + <)-^ 

with bniO = \/2sin(7rO- The eigen-pairs (vp, (f)p{^)) are defined by the following 
differential equation and the boundary condition: 



d2(^„ / o 



d^2 I (1 + a^)2 



, <Pp{0) = <j)p{l) = . (D.7) 



We can easily recognize that the solutions of equation (D.7) converge in the limit a ^ 
to Vp = hp and 0^(0 = bp{^). We assume that the solutions can be expanded in a 
power series of variable a. The eigen-pairs can then be obtained using the standard 
perturbation theory with the perturbation parameter a. The rescaled mode numbers 
are written as 



2 u2 

= hp 



1 + a + aU - 



1 1 



a'(\ + 0{hi))+0{a')^ (D. 



.3 2(7rj9)2^ 
and the rescaled mode functions read as 

Snn hj^ 

m) = Y.^nMO, = ^np + « ^4(X ^2)3 Ct + 0{o?) , (D.9) 
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where we use the symbol 5^^^ = (1 : n is odd; : otherwise). By plugging the mode 
functions (j)p{C) (D.8) into transition matrices Q (D.5) and P (D.6) we obtain 



Note that the rescaled transition matrices Q and P satisfy the known identity QP^ = 
PQ^ = id. We insert the expressions for Q (D.IO) and P (D.ll) back into S± (D.3) and 
C± (D.2) and write the reflection and the transmission matrix as 



where we have introduced the symbol (A)^ = A + A^ and the diagonal matrix 
H = didig{hn}n=i- The approximations of the reflection matrix R (D.13) and the 
transmission matrix T (D.14) are valid far away from the resonant condition ak = 7rn, 
because we assumed that l^'^l > a for all m < Nq. We conclude that the strength 
of reflection scales as 7^ ^ and that narrow channels can be treated as perturbed 
straight wave-guides. 




(D.IO) 



(D.ll) 



(D.12) 




